{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "4fbae0f9",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "# Molecule Eigensolver (VQE method)\n",
    "\n",
    "Evaluating the ground state of a molecular Hamiltonian allows us to understand the chemical properties of the molecule. In this demo, we demonstrate the usage of Variational Quantum Eigensolver (VQE) for finding the ground states and energies of several molecules:  \ud835\udc3b2 ,  \ud835\udc3b2\ud835\udc42  and  \ud835\udc3f\ud835\udc56\ud835\udc3b .\n",
    "\n",
    "VQE is a leading method for finding approximate values of ground state wavefuncions and energies for complicated quantum systems, and as such can give solutions for complex molecular structures. The overview of the VQE method is as following: a problem (i.e. a molecule) is defined by a Hamiltonian which ground state is sought. Then, a choice of a parameterized ansatz is made. Using a hybrid quantum-classical algorithm, a solution for the defined parameters which minimize the expectation value for the energy is found. A clever ansatz will lead to an estimated ground state solution.\n",
    "\n",
    "Within the scope of Classiq's VQE algorithm, the user defines a Molecule which is being translated to a concise Hamiltonian. Then, a choice between several types of well studied ansatz is given, which can be carefully picked to fit your Molecule type. In the last stage, the Hamiltonian and ansatz are sent to a classical optimizer. During this tutorial we will demonstrate the steps and user options in Classiq's VQE algorithm. Furthermore, the demonstration will present the optimization strength Classiq's VQE algorithm holds, and it's state of the art results in term of efficient quantum circuit - with ultimate combination of low depth and high accuracy, while minimizing the number of CX gates."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "54d09062-1b3b-4e4b-8351-5e05a633e269",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "## 0. Pre-requirments\n",
    "\n",
    "The model is using several Classiq's libraries."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "c6bbe65f-1e6a-475c-a43f-cb4cc04bbdfa",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2024-05-07T15:21:54.983816Z",
     "iopub.status.busy": "2024-05-07T15:21:54.983159Z",
     "iopub.status.idle": "2024-05-07T15:21:58.305203Z",
     "shell.execute_reply": "2024-05-07T15:21:58.304507Z"
    },
    "pycharm": {
     "name": "#%%\n"
    }
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "\n",
    "from classiq import QuantumProgram, construct_chemistry_model, execute, show, synthesize\n",
    "from classiq.applications.chemistry import (\n",
    "    ChemistryExecutionParameters,\n",
    "    HEAParameters,\n",
    "    Molecule,\n",
    "    MoleculeProblem,\n",
    "    UCCParameters,\n",
    ")\n",
    "from classiq.execution import (\n",
    "    ClassiqBackendPreferences,\n",
    "    ClassiqSimulatorBackendNames,\n",
    "    ExecutionPreferences,\n",
    "    OptimizerType,\n",
    ")\n",
    "from classiq.synthesis import set_execution_preferences"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "faa3c10f",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "## 1. Generate Qubit Hamiltonian\n",
    "\n",
    "The first step is to define the molecule we wish to simulate. We hereby declare the class Molecule and insert a list of atoms and their spacial positions. The algorithm will automatically regard relevant attributes such as the atom's mass, charge and spin.\n",
    "\n",
    "As mentioned above, during this tutorial, we demonstrate how to define and find the ground state and energies for 3 molecules:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "ef2ce57e",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2024-05-07T15:21:58.308346Z",
     "iopub.status.busy": "2024-05-07T15:21:58.307789Z",
     "iopub.status.idle": "2024-05-07T15:21:58.312947Z",
     "shell.execute_reply": "2024-05-07T15:21:58.312497Z"
    },
    "pycharm": {
     "name": "#%%\n"
    }
   },
   "outputs": [],
   "source": [
    "molecule_H2 = Molecule(atoms=[(\"H\", (0.0, 0.0, 0)), (\"H\", (0.0, 0.0, 0.735))])\n",
    "molecule_O2 = Molecule(atoms=[(\"O\", (0.0, 0.0, 0)), (\"O\", (0.0, 0.0, 1.16))])\n",
    "molecule_LiH = Molecule(atoms=[(\"H\", (0.0, 0.0, 0.0)), (\"Li\", (0.0, 0.0, 1.596))])\n",
    "molecule_H2O = Molecule(\n",
    "    atoms=[(\"O\", (0.0, 0.0, 0.0)), (\"H\", (0, 0.586, 0.757)), (\"H\", (0, 0.586, -0.757))]\n",
    ")\n",
    "molecule_BeH2 = Molecule(\n",
    "    atoms=[(\"Be\", (0.0, 0.0, 0.0)), (\"H\", (0, 0, 1.334)), (\"H\", (0, 0, -1.334))]\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ab162f48",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "Similarly, the user is able to construct any valid essambly of atoms. The distances are recived in \u00c5 ($10^{-10} m$). We will continue this demonstration with a specific molecule. The user can change the `molecule` below to study other cases."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "2e77678f",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2024-05-07T15:21:58.315498Z",
     "iopub.status.busy": "2024-05-07T15:21:58.314960Z",
     "iopub.status.idle": "2024-05-07T15:21:58.318124Z",
     "shell.execute_reply": "2024-05-07T15:21:58.317557Z"
    },
    "pycharm": {
     "name": "#%%\n"
    }
   },
   "outputs": [],
   "source": [
    "molecule = molecule_H2"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c4a541a5",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "Next, we define the parameters of the Hamiltonian generation program. The user has a choice over the following options:\n",
    "- mapping (str): the mapping between the fermionic Hamiltonian and an qubits Hamiltonian. Supported types:\n",
    "        - \"jordan_wigner\"\n",
    "        - \"parity\"\n",
    "        - \"bravyi_kitaev\"\n",
    "        - \"fast_bravyi_kitaev\"\n",
    "- freeze_core (bool): remove the \"core\" orbitals of the atoms defining the molecule.\n",
    "- z2_symmetries (bool): whether to perform z2 symmetries reduction. If symmetries in the molecules exist, this option will decrease the number of qubits used and will efficient the Hamiltonian and thus the calculations.\n",
    "\n",
    "Finally, the Hamiltonian is generated from `MoleculeProblem`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "2e0426d5",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2024-05-07T15:21:58.320678Z",
     "iopub.status.busy": "2024-05-07T15:21:58.320142Z",
     "iopub.status.idle": "2024-05-07T15:22:03.633083Z",
     "shell.execute_reply": "2024-05-07T15:22:03.632220Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Your Hamiltonian is\n",
      "-1.041 * I\n",
      "-0.796 * Z\n",
      "+0.181 * X\n"
     ]
    }
   ],
   "source": [
    "chemistry_problem = MoleculeProblem(\n",
    "    molecule=molecule,\n",
    "    mapping=\"jordan_wigner\",  #'bravyi_kitaev'\n",
    "    z2_symmetries=True,\n",
    "    freeze_core=True,\n",
    ")\n",
    "\n",
    "operator = chemistry_problem.generate_hamiltonian()\n",
    "gs_problem = chemistry_problem.update_problem(operator.num_qubits)\n",
    "print(\"Your Hamiltonian is\", operator.show(), sep=\"\\n\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "67a124c2-33fc-420c-8cb3-0002f18f4e77",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "The output of the above code lines is the Hamiltonian presented as a superposition of Pauli matrices multiplication.\n",
    "One can easily confirm that using z2*symmetries=True, the number of qubits are reduced (compered to z2_symmetries=False): for $H_2$ - from 4 to 1, for $LiH$ from 12 to 8, and for $H*{2}O$ from 14 to 10."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8abe3e3d-1b01-4fab-b86a-feaab3851950",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "## 2. Constructing and Synthesizing a Ground State Solver\n",
    "\n",
    "A ground state solver model consists of a parameterized eigenfunction (\"the ansatz\"), on which we run a VQE. In addition, a post-process of the result allows to return the total energy (combining the ground state energy of the Hamiltonian, the nuclear repulsion and the static nuclear energy).\n",
    "\n",
    "Once we've specified an Hamiltonian and a desired Ansatz, we send them to the VQE algorithm in order to find the Hamiltonian's ground state. In the process, the algorithm will send requests to a classical server, which task is to minimize the energy expectation value and return the optimized parameters. The simulator and optimizing parameters are defined as part of the VQE part of the model. The user should control the `max_iteration` value in a manner so the solution has reached a stable convergence. In addition, the value `num_shots` sets the number of measurements performed after each iteration, thus influence the accuracy of the solutions.\n",
    "\n",
    "We demonstrate two different proposal for the wavefunction solution ansatz: (1) Hardware (HW) efficient, and (2) Unitary Coupled Cluster (UCC). For groundstate solvers it is typical to initialize the Ansatz with the Hartree-Fock state."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "07ed8f8e-485b-4f0f-93b7-78e7fd6fbe39",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "### 2.1 HW-Efficient Ansatz\n",
    "\n",
    "Hardware-efficient ansatz is a suggested solution that is generated to fit a specific hardware [1]. The ansatz creates a state with given number of parameters by user choice (number of qubits, that should fit the Hamiltonian), and creates entanglement between the qubits by the inputed connectivity map. In this example, a 4 qubit map is given, which is specifically made of $H_2$ with z2_symmetries=False.\n",
    "\n",
    "After constructing the model, we can synthesize it and view the outputted circuit, in charged on creating the state with an interactive interface."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "90b20061-8dbd-4136-adba-28ddacb1f583",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2024-05-07T15:22:03.637345Z",
     "iopub.status.busy": "2024-05-07T15:22:03.637033Z",
     "iopub.status.idle": "2024-05-07T15:22:03.660043Z",
     "shell.execute_reply": "2024-05-07T15:22:03.659362Z"
    }
   },
   "outputs": [],
   "source": [
    "chemistry_problem = MoleculeProblem(\n",
    "    molecule=molecule,\n",
    "    mapping=\"jordan_wigner\",  #'bravyi_kitaev'\n",
    "    z2_symmetries=False,\n",
    "    freeze_core=True,\n",
    ")\n",
    "\n",
    "hwea_params = HEAParameters(\n",
    "    num_qubits=4,\n",
    "    connectivity_map=[(0, 1), (1, 2), (2, 3)],\n",
    "    reps=3,\n",
    "    one_qubit_gates=[\"x\", \"ry\"],\n",
    "    two_qubit_gates=[\"cx\"],\n",
    ")\n",
    "\n",
    "qmod = construct_chemistry_model(\n",
    "    chemistry_problem=chemistry_problem,\n",
    "    use_hartree_fock=True,\n",
    "    ansatz_parameters=hwea_params,\n",
    "    execution_parameters=ChemistryExecutionParameters(\n",
    "        optimizer=OptimizerType.COBYLA,\n",
    "        max_iteration=30,\n",
    "        initial_point=None,\n",
    "    ),\n",
    ")\n",
    "\n",
    "backend_preferences = ClassiqBackendPreferences(\n",
    "    backend_name=ClassiqSimulatorBackendNames.SIMULATOR\n",
    ")\n",
    "\n",
    "qmod = set_execution_preferences(\n",
    "    qmod,\n",
    "    execution_preferences=ExecutionPreferences(\n",
    "        num_shots=1000, backend_preferences=backend_preferences\n",
    "    ),\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "f13599d2-e8e7-4165-aa26-8ae8d4bdafaa",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2024-05-07T15:22:03.663419Z",
     "iopub.status.busy": "2024-05-07T15:22:03.663167Z",
     "iopub.status.idle": "2024-05-07T15:22:03.686164Z",
     "shell.execute_reply": "2024-05-07T15:22:03.685472Z"
    }
   },
   "outputs": [],
   "source": [
    "from classiq import write_qmod\n",
    "\n",
    "write_qmod(qmod, name=\"molecule_eigensolver\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "22cd12d1-2c87-400a-a983-b2f24e40fa45",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2024-05-07T15:22:03.690648Z",
     "iopub.status.busy": "2024-05-07T15:22:03.689438Z",
     "iopub.status.idle": "2024-05-07T15:22:17.207052Z",
     "shell.execute_reply": "2024-05-07T15:22:17.206288Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Opening: https://platform.classiq.io/circuit/dd765d79-bcdb-49e0-85be-74af8c1b3f86?version=0.41.0.dev39%2B79c8fd0855\n"
     ]
    }
   ],
   "source": [
    "qprog = synthesize(qmod)\n",
    "show(qprog)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7fb92cdc",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "### 2.2. UCC Ansatz\n",
    "\n",
    "Next, we show how to create the commonly used chemistry-inspired UCC ansatz, which is a unitary version of the classical coupled cluster (CC) method [2] .\n",
    "\n",
    "The parameter that defines the UCC ansatz is:\n",
    "- excitations (List[int] or List[str]): list of desired excitations. Allowed excitations:\n",
    "        - 1 for singles\n",
    "        - 2 for doubles\n",
    "        - 3 for triples\n",
    "        - 4 for quadruples\n",
    "\n",
    "Once again, after the code lines bellow run, the user is able to view the outputted circuit, in charged on creating the state with an interactive interface. In addition, the depth of the circuit is printed.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "1f520673",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2024-05-07T15:22:17.210966Z",
     "iopub.status.busy": "2024-05-07T15:22:17.210537Z",
     "iopub.status.idle": "2024-05-07T15:22:42.664650Z",
     "shell.execute_reply": "2024-05-07T15:22:42.663927Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Opening: https://platform.classiq.io/circuit/39cf374d-1965-4081-85bf-75c71cd73296?version=0.41.0.dev39%2B79c8fd0855\n",
      "circuit depth: 3\n"
     ]
    }
   ],
   "source": [
    "chemistry_problem = MoleculeProblem(\n",
    "    molecule=molecule,\n",
    "    mapping=\"jordan_wigner\",  #'bravyi_kitaev'\n",
    "    z2_symmetries=True,\n",
    "    freeze_core=True,\n",
    ")\n",
    "\n",
    "serialized_chemistry_model = construct_chemistry_model(\n",
    "    chemistry_problem=chemistry_problem,\n",
    "    use_hartree_fock=True,\n",
    "    ansatz_parameters=UCCParameters(excitations=[1, 2]),\n",
    "    execution_parameters=ChemistryExecutionParameters(\n",
    "        optimizer=OptimizerType.COBYLA,\n",
    "        max_iteration=30,\n",
    "        initial_point=None,\n",
    "    ),\n",
    ")\n",
    "\n",
    "backend_preferences = ClassiqBackendPreferences(\n",
    "    backend_name=ClassiqSimulatorBackendNames.SIMULATOR\n",
    ")\n",
    "\n",
    "serialized_chemistry_model = set_execution_preferences(\n",
    "    serialized_chemistry_model,\n",
    "    execution_preferences=ExecutionPreferences(\n",
    "        num_shots=1000, backend_preferences=backend_preferences\n",
    "    ),\n",
    ")\n",
    "\n",
    "qprog = synthesize(serialized_chemistry_model)\n",
    "show(qprog)\n",
    "\n",
    "circuit = QuantumProgram.from_qprog(qprog)\n",
    "print(f\"circuit depth: {circuit.transpiled_circuit.depth}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f12c9a32-f271-4892-8515-bd9e6b9fcf8b",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "Classiq's UCC algorithm provides an highly efficient solution in aspects of circuit depth and number of CX gates. Those ultimately reduce the gate's time and amount of resources needed for its operation."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ef36661f",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "## 3. Execute to Find Ground State\n",
    "\n",
    "Once we've synthesized the model we can execute it."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "1a66d377",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2024-05-07T15:22:42.670417Z",
     "iopub.status.busy": "2024-05-07T15:22:42.669946Z",
     "iopub.status.idle": "2024-05-07T15:22:50.224366Z",
     "shell.execute_reply": "2024-05-07T15:22:50.223567Z"
    }
   },
   "outputs": [],
   "source": [
    "result = execute(qprog).result()\n",
    "chemistry_result_dict = result[1].value"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f72ceeb5-bc71-46d8-a390-31ec874700f3",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "Execution of the quantum program returns several useful outputs:\n",
    "- energy : the output of the VQE algorithm - the electronic energy simulated.\n",
    "- nuclear_repulsion : the electrostatic energy generated by the atom's nuclei.\n",
    "- hartree_fock_energy : the Hartree Fock energy.\n",
    "- total_energy : this is the ground state energy of the Hamiltonian (combining the energy, the nuclear repulsion and the static nuclear energy).\n",
    "\n",
    "It also contains the full VQE result from which we can get, for example:\n",
    "- optimal_parameters : gives the results for the anzats parameters minimizing that expectation value.\n",
    "- eigenstate : gives the ground state wave function.\n",
    "\n",
    "Note the all energy are presented in units of Hartree."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "437b3211",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2024-05-07T15:22:50.229659Z",
     "iopub.status.busy": "2024-05-07T15:22:50.228329Z",
     "iopub.status.idle": "2024-05-07T15:22:50.239525Z",
     "shell.execute_reply": "2024-05-07T15:22:50.238826Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "-1.1382090268147285"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "chemistry_result_dict[\"total_energy\"]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "id": "9a537d3c",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2024-05-07T15:22:50.244086Z",
     "iopub.status.busy": "2024-05-07T15:22:50.242988Z",
     "iopub.status.idle": "2024-05-07T15:22:50.250371Z",
     "shell.execute_reply": "2024-05-07T15:22:50.249692Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "{'param_0': 3.2624618497328566}"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "chemistry_result_dict[\"vqe_result\"][\"optimal_parameters\"]"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2375f3c3",
   "metadata": {},
   "source": [
    "Finally, we can compare the VQE solution to the classical solution by employing exact diagonalization:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "id": "5c896576",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2024-05-07T15:22:50.255196Z",
     "iopub.status.busy": "2024-05-07T15:22:50.253856Z",
     "iopub.status.idle": "2024-05-07T15:22:50.261974Z",
     "shell.execute_reply": "2024-05-07T15:22:50.261293Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "exact result: -1.8572750302023786\n",
      "vqe result: -1.8581780212637082\n"
     ]
    }
   ],
   "source": [
    "mat = operator.to_matrix()\n",
    "w, v = np.linalg.eig(mat)\n",
    "print(\"exact result:\", np.real(min(w)))\n",
    "print(\"vqe result:\", chemistry_result_dict[\"energy\"])"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b998e321",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "[1] Abhinav Kandala, Antonio Mezzacapo, Kristan Temme, Maika Takita, Markus Brink, Jerry M. Chow, Jay M. Gambetta Hardware-efficient variational quantum eigensolver for small molecules and quantum magnets. Nature 549, 242 (2017)\n",
    "\n",
    "[2] Panagiotis Kl. Barkoutsos, Jerome F. Gonthier, Igor Sokolov, Nikolaj Moll, Gian Salis, Andreas Fuhrer, Marc Ganzhorn, Daniel J. Egger, Matthias Troyer, Antonio Mezzacapo, Stefan Filipp, and Ivano Tavernelli Quantum algorithms for electronic structure calculations: Particle-hole Hamiltonian and optimized wave-function expansions Phys. Rev. A 98, 022322 (2018)\n"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.11.9"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}
